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<N i ABSTRACT 



A magnetic flux tube may be considered both as a separate body and as a 
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. confined field. As a field, it is affected both by differential rotation (f2-effect) and 

cyclonic convection (cc-effect) . As a body, the tube experiences not only a buoy- 
ant force, but also a dynamic pressure due to downflows above the tube. These 
two competing dynamic effects are incorporated into the a-Q dynamo equations 
through the total magnetic turbulent diffusivity, leading to a flux tube dynamo 
operating in the convection zone. We analyze and solve the extended dynamo 
equations in the linear approximation by adopting the observed solar internal ro- 
tation and assuming a downflow effect derived from numerical simulations of solar 
convection zone. The model reproduces: the 22-year cycle period; the extended 
butterfly diagram with the confinement of strong activity to low heliographic 
latitudes |$| < 35°; the evidence that at low latitudes the radial field is in an 
approximately 7r phase lag compared to the toroidal field at the same latitude; 
the evidence that the poleward branch is in a tt/2 phase lag with respect to the 
equatorward branch; and the evidence that most of the magnetic flux is present 
in an intermittent form, concentrated into strong flux tubes. 

Subject headings: Sun: convection zone — Sun: magnetic fields — Sun: dynamo 
theory 



1. Introduction 



The magnetic activity of the Sun and late type stars is currently attributed to a dynamo 
action, namely the interaction of conducting fluid motions and magnetic fields in a suitable 
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spherical shell. This was first suggested by Larmor (1919), but the development of these ideas 
proved to be very complex, especially in view of Cowling's theorem (1934), which indicated 
that a dynamo could not operate with strictly axisymmetric magnetic fields. The generation 
and evolution of magnetic fields in the Sun is mostly described in terms of an a-Q dynamo 
operating in the convection zone or in the boundary (overshoot) layer just beneath the 
convection zone, as the well known recent results of helioseismology may suggest. Although 
the physics of the interaction of rotation and convective modes is very complex, and poorly 
understood in spite of ever more refined attempts, the simple (in principle) approach based on 
differential rotation and cyclonic turbulence has been productive in modeling solar activity. 

A good dynamo model should, of course, reproduce the main characteristics of the 
solar cycle: equatorward migration of spots and, in general, of active regions where strong 
magnetic fields are observed, starting from a latitude of about 35°; poleward migration of 
some magnetic tracers such as high latitude prominences and polar faculae (e.g., Makarov 
& Sivaraman 1989); polarity reversals of the polar fields every 11-years; polarity rules in the 
two hemispheres (Hale's law); ratio and phase relation of the poloidal field to the toroidal 
field ((Stix 1976)); 22- year cycle period; intermittent form of most of the magnetic flux, 
which is mainly concentrated into flux tubes; rigid rotation of the preferential longitudes 
(some fixed longitudes where spots show a rigid rotation); flip-flop effect in the preferential 
longitudes (namely the alternate occurrence of spots at these longitudes in the two 11-year 
halves of the 22-year magnetic cycle); and possibly long-term periodicities shown by solar 
activity, and absence of activity at some epochs (Maunder minima). 

The location of the a-Q dynamo action is an important theoretical problem. In the 
course of the years, several locations have been proposed: mostly the whole convection 
zone or the stable layers at the bottom of the convection zone (this has been done in 
an extensive series of papers starting from the fundamental work by Parker (1955)); but 
also the boundary (overshoot layer) just beneath the convection zone (e.g., Belvedere et 
al. 1991, Rudiger & Brandenburg 1995). It is usually assumed that the and a-effects 
are located in the same region, since a spatial separation is not plausible as would im- 
ply problems of upward and downward magnetic field transport which are not easily over- 
come by means of the turbulent diffusion. An exception are the so-called interface mod- 
els (Parker 1993; Charbonneau & MacGregor 1997; Markiel & Thomas 1999), that assume 
the f2-effect to operate on the low diffusivity side of an interface (the base of the con- 
vection zone), while the a-effect works on the high diffusivity side, and, in some aspects, 
also the advective models (Choudhuri et al. 1995; Durney 1995; Durney 1996; Durney 1997; 
Nandy & Choudhuri 2001), see also the recent work by Bonanno et al. (2002). In these mod- 
els, the location of the two effects may be varied, sometimes having the characteristics of an 
interface dynamo, sometimes of a Babcock-Leighton dynamo (with a deeply located f2-effect 
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and a surface located ct-effect). The advection dynamos essentially work like a conveyor 
belt that is provided by the meridian circulation in the form of a single or more cells in the 
convection zone (Wang & Sheeley 1991; Dikpati & Charbonneau 1999). 

Of course the validity of the various approaches and models is measured by their ca- 
pability to reproduce as many features of the solar activity as possible. But this is not the 
only criterion to evaluate a dynamo model: in fact, the current framework of solar dynamo 
theory, which appears to be in agreement with the basic features of solar activity, arose from 
a series of attempts that dealt separately with partial, although fundamental aspects of the 
problem. 

Among those attempts, in the framework of this paper, we consider flux tube dynamos. 
It is likely that, in future work on MHD of stellar interiors and dynamo theory, flux tube 
dynamics will be particularly important. Earlier studies on the dynamics of thin flux tubes, 
intermittently erupting at the solar surface, were carried out by (e.g.) Schuessler (1980), 
Spruit (1981), Spruit and Van Ballegooijen (1982a,b), van Ballegooijen (1983), Schmitt 
(1987), Ferriz-Mas and Schuessler(1994,1995), Caligari et al.(1995), and Schuessler et al. 
(1996). 

In the thin tube approximation, the tube radius must be much smaller than its radius 
of curvature, the local scale height, the rotational shear scale, and the size of the relevant 
perturbations, including the large scale turbulence size (i.e. the size of dominant eddies). 
The basic point is that strong axisymmetric toroidal fields located in the subadiabatic and 
relatively stable layers at the bottom of the convection zone, generated from poloidal fields 
by rotational shear, can be represented as an ensemble of isolated flux tubes, subjected to 
magnetic buoyancy and other perturbations like aerodynamical drag, downward flows and 
rotationally induced effects such as Coriolis force and differential centrifugal force. 

In this paper, we argue that, due to the competition between two effects, namely 
buoyancy force and downward motions, the rising time of the tubes may be comparable 
to the timescale of the solar cycle. In fact, numerical simulations of the outer solar con- 
vection zone indicate a major presence of downward moving plumes (Chan and Sofia 1989; 
Kim & Chan 1998; Nordlund 1999; Robinson et al 2003). A similar idea was already pre- 
sented by Parker (1975), who considered the aerodynamical drag as the competing effect. 
The new approach we adopt in this paper consists in giving a particular relevance to the 
downward flows and in incorporating the two competing effects into the mean field dynamo 
equations, by deriving an expression for the total turbulent diffusivity that takes into ac- 
count the action of both buoyancy and downward motions. Then, by solving the linear 
problem in the usual perturbative methods, we obtain, in a self-consisting way, results in 
close agreement with the main observed solar cycle features. A good feature of this model is 
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the self-consistency of the approach formulated in the framework of a single mathematical 
(and physical) development, without any ad hoc parameterizations. 

In this model the solar internal differential rotation deduced from inversions of helio- 
seismic observations (e.g., Antia et al. 1998) has been adopted. 

We show in §2 how the effects of magnetic buoyancy and downward motions can be 
incorporated into the mean field dynamo equations by modifying the usual expression for 
the turbulent diffusivity (5. In §3 we linearly analyze the dynamo equations to obtain the 
observable quantities of the model. In §4 we use the theoretical findings to interpret the 
properties of the solar cycle listed above. Finally, in §5 we describe the main conclusions of 
the paper. 

Notice that in this analysis antisymmetrical (dipole-like) modes are naturally assumed, 
as observations show that they are dominant in the Sun. Indeed, it is beyond the scope 
of this paper to solve a bidimensional (r, 9) eigenvalue problem in order to check whether 
dipole- modes are more easily excited and why (as was done, e.g., by Bonanno et al. 2002). 



2. Basic equations 

We consider a magnetic flux tube with radius a and length L. Its magnetic field B is 
assumed to be uniform. As usual, we decompose it into a toroidal and a poloidal component, 

B = Be <f) + V x (Ae ), (1) 

where is the azimuthal unit vector, B is the toroidal field component, and A is the 
(toroidal) vector potential of the poloidal field. The local hydrostatic equilibrium of the flux 
tube in the solar interior requires 

P = P t + B 2 /8n, (2) 

where P and Pj are the total and internal gas pressure, respectively. Assuming that this 
equilibrium in the flux tube is reached solely by adjusting its density, then the internal 
density pi is related to the normal density p by 

Pi = pPi/P. (3) 

As a result, the density reduction of the magnetic flux tube is 



p - p t = P B 2 /8nP. 



(4) 
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The reduced density inside the flux tube produces a buoyancy, 

F = na 2 Lg(p - pi)e r , (5) 

where g is the gravitational acceleration. We assume that a is smaller than the pressure 
scale height so that we can consider that the gas density of the flux tube is uniform. In this 
case, the total mass contained in the flux tube equals 

Mi = ira 2 Lp(l-B 2 /B g ), (6) 

where B g = (8-kP) 1 ^ 2 . Therefore, the buoyant acceleration g b can be expressed by the 
gravitational acceleration g, total pressure P, and the magnetic field strength B, as follows: 

g b = gB 2 /(B 2 -B 2 ). (7) 

If there is no downward flow to balance this pressure, the only force that goes against 
the buoyant force is the aerodynamic drag, 

F d = -\C d pu 2 aLe r , (8) 

where u is the velocity of rise and Cd(~ 1) is the drag coefficient. Taking into account both 
the buoyant force and the drag force, as done by Parker (1975), the terminal velocity of rise 
occurs for F = F D , yielding the rate of rise 

o \l/2 / \ 1/2 

Zirag p — Pi\ I rra 



where H p = P/pg is the pressure scale height, Va = B/(Anp) 1 ^ 2 is the Alfven speed. The 
rise time to the surface is smaller than the required field amplification time of about 5 years 
by the classical dynamo theory, as estimated by Parker (1975). 

However, numerical simulations of the Sun's outer convection zone made, e.g. by Chan 
& Sofia (1989), indicate a major presence of downward-moving plumes of high velocity. This 
downward flow has an inward dynamic pressure pv 2 . If this pressure can balance the buoyant 
force of the flux tube, then the flux tube can be stored in the convection zone so that it 
can be amplified by the a-fl dynamo. As a body, a magnetic flux tube experiences not only 
a buoyancy, but also a dynamic pressure of the downflows. The total pushdown force of 
downward plumes equals to the dynamic pressure {pv 2 ) times the total cross section of the 
plumes (S). The resulting acceleration g v can be expressed as follows: 



g v = 2sv 2 M™(l-B 2 /B 2 )}, 



(10) 
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where s — 1 when v z > 0, s = — 1 when < 0, £ = S/2aL < 1 is the fractional area of 
the downflows. The upward-moving plumes underneath the flux tube (s = 1) accelerate the 
magnetic buoyant diffusion of the tube, while the downward-moving plumes above the flux 
tube (s = — f ) go against the magnetic diffusion. These downflows act as an anti-diffusion 
process. In general, — 1 < s < +1, depending on the localtion of the flux tube. 

The condition that the dynamic pressure of the downflows balances the buoyancy de- 
termines a critical value of the magnetic field B c , which is given by + g v = 0, 

B c = B g {-2svlil-Kagfl\ (II) 

When B < B c , the buoyancy is smaller than the dynamic pressure of the downflows, the 
dynamo magnetic field may grow; when B > B c , the dynamo field may damp. These two 
ingredients, buoyancy and downflows, are shown in Fig. 1. Since downflows occur in the 
convection zone and the overshoot region near the base of the convection zone, this dynamo 
mechanism may operate in both of these regions. 

The magnetic field in a flux tube can be considered to be a "mean field" since the 
mean-field concept is relative (e.g., Krause and Radler 1980). In fact, it depends on the 
spatio-temporal range over which the average is made. For the mean field in a flux tube, the 
spatio-temporal range is the volume and lifetime of a typical magnetic flux tube. With this 
in mind, the field in a flux tube obeys the classic dynamo equation (Parker 1955), 

- = Vx(UxB + ^)-Vx(i)VxB). (12) 

Using this equation, we cannot study the formation of flux tubes from the large-scale mag- 
netic field, but we can determine if the flux tube fields can be amplified. We assume a pure 
rotation motion U = Q x r, where f2 is the angular velocity. This will cause the Q-effect. 

Both the buoyant and dynamic diffusivities are direction-dependent. So they are anisotropic. 
Such anisotropic diffusivities must be treated as tensors in principle (Krause & Radler 1980). 
However, it is hard to do so in practice. To our knowledge, almost all dynamo models (in- 
cluding most conveyor belt ones) consider thermal, dynamical and magnetic diffusivities as 
scalars. All recent books that deal also with solar dynamo models (e.g., Stix 1989; Mestel 
1999; Ruediger & Hollerbach 2004) quote and review only papers that assume scalar dif- 
fusivities. So, a solar model that includes the anisotropy of diffusivity in tensorial form is 
currently beyond the state of the art. 

Physically, the tensorial form of diffusivity should be used only when a strong source 
of asymmetriy exists. This is not the stellar case, where spherical symmetry is generally 
assumed, and no strong axial electromotive force, no fast rotation occur. Practically, the 
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dynamo theory that uses the scalar form of cliff usivities has captured the essentials of the 
solar cycle, and is able to reproduce the butterfly diagram, and other observational facts. 
Mathematically, the scalar form of diffusivities can be considered to be the zero-order approx- 
imation of the general tensorial form. For these reasons, we use the following approximation 
form of the general electromotive force (Krause & Radler 1980) 

^-«( 0) B-/3 (0) VxB, (13) 

where 

/3 (0) KV + Pt + Pb + Pv (14) 

Defining 

a = -a {0 \ (3 = (3 {0) , (15) 
we obtain the governing equation of a magnetic flux tube, 

= V x (U x B + aB) - V x (/5V x B). (16) 

at 

Since (5 depends on B, this equation is nonlinear. In the azimuthal symmetric case, this 
equation is equivalent to the following two scalar equations: 

dA „{ldA cotOdA d 2 A Id 2 A 

aB + (3( + 



dt \r dr r 2 89 dr 2 r 2 d6 2 

(17) 

OB A dA dA (d 2 A 1 d 2 A s 

a x A + a 2 ^ — V «3tt7t _ a ^rr + 



dt dr ° 89 \dr 2 r 2 dd 2 



, „ dB dB a (d 2 B 1 d 2 B 



(d 2 B 1 d 2 B\ 

[Ik 2 + ^ 2 ~dl 2 ) ' 

(18) 



See appendix A for the detailed derivation and the definitions for those coefficients a s and 
b s . Here B and A are defined in Eq. (1). 

The magnetic diffusion caused by the microscopic gas collision motion, has been included 
in Eq. (12) in terms of the classic magnetic diffusivity 77 V x B. For a homogeneous convection, 
a can be approximately expressed as (Krause & Radler 1980) 

ol = \hr COT , (19) 



where h = v • (V x v) is the mean helicity over the correlation time r cor . But in the general 
case, a may deviate from this expression. Since the turbulent velocity v is related to the 



- 8- 



background and dynamo magnetic field B, a depends upon B. This means that a-effects 
must be nonlinear. 

In order to work out the various magnetic diffusivities mentioned above, we observe 
that the diffusivity of a force is proportional to its impulse per unit mass. The impulse is 
defined as the force times its action time. The impulse per unit mass has the dimension 
of acceleration times time. The action time can be taken to be the dynamic response time 
of a magnetic tube. Since it is assumed that the process is iso-thermal in calculating the 
buoyant acceleration, the dynamic response time is the travel time of a perturbation across 
the diameter of the flux tube at the isothermal sound speed, C s = (P/p) 1//2 , so t act = 2a/C s . 
However, diffusivity has the dimension of acceleration times time times length. The shortest 
characteristic length for a flux tube is its diameter a. We use this length a to meet the need 
of dimensionality for the buoyant and the dynamic diffusivity. Consequently, we have 

Pb = #/^ ^ , (20) 

P v = g v (2a/C s )2a = ^{B 2 - B 2 ) ^ ' 

The magnetic turbulent diffusivity (3 t = \v ' 2 r cor for an isotropic turbulence is known (Krause & Radler 1980 
where v' is the turbulent velocity. 



3. Linear analysis 

The flux tube model equations obtained in the previous section are nonlinear partial- 
differential equations. Since such equations are difficult to solve analytically, while a simple 
linear analysis may reveal some important characteristics about their solution, we analyze 
and solve the equation in the linear approximation in this paper. Unlike numerical simula- 
tions, the linear analysis cannot take into account full nonlinearity though it reveals some 
nonlinear effects such as critical magnetic fields of a flux tube. 

For the sake of simplicity, we neglect rj in this analysis, namely, we assume 

P^Pt + Pb + Pv (22) 
The second simplification assumption is B 2 <C B g , which means 



Pb « ( P /P) 1/2 4a 2 gB 2 /B 2 g , 
P v » (p/PY^asv^/n. 



(23) 
(24) 
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3.1. Linearized dynamo equations 

Since b s and (3 in Eqs. (17-18) depend upon B, dB/dr and/or dB/89, we need to 
specify B = Bo and A = Aq as the reference state to investigate their perturbations B' 
and A' near the reference state. See Appendix B for the detailed linearization of Eqs. (17- 
18). Omitting all the superscripts and subscripts , and bearing in mind that all the 
coefficients are evaluated at the reference state, we obtain the linearized dynamo equations 
from Eqs. (B3-B4) when the a 2 -effects are neglected: 



where 



OA' 
~W 



= (a + c )B' 

| | cotedA' | d 2 A' | 1 d 2 A'\ 

y r Q r r 2 QQ Q r 2 r 2 QQ2 J ' 

8B' Al dA' dA 

~dt = aiA+a2 ^ + a3 W 

OB' OB' n fd 2 B' 1 8 2 B' 
+dB' + c 2 — + c 3 — - + (3 -— + 



dr 86 V dr 2 r 2 d9 



d\ = — (cos 9& r — Cl e sin 9), 

r 

a 2 = —QCIq sin 9, 

Q * 

a 3 = — f2 r sin 9, 

r 

2f% (IdA cot9_dA cPA 1 d 2 A\ 
C ° ~ ~B \r~dr~ + ~^ 2 ~~d9 + W 2 + 1*~W ) ' 

ci = ^(AB r + 3B e cot9)-^(2(3 b + (3)csc 2 9 + ^ 
2A,a 2 , a 2 ^ , 2(3 b fd 2 B , 1 d 2 B 



(26) 



ry ry ry 

c 3 = ^cot0 + ^(cot0 + 2B fl ). 

We define f2 r , fig, and (3 is defined in Appendix A. Similarly, B r = dlnB/dhir, Bg = 
d\nB/d9. 
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3.2. Dispersion relation 

In order to investigate the linear instability, we try the wavelike solutions of the form 

B' = B£exp[i(wf-k-x)], (27) 
A' = A' exp[i{ut-k-x)}, (28) 

where u> is complex, while the other variables are real. Substituting Eqs. (27-28) into Eqs. (25- 
26), we obtain 

[iuj + (3(ik r /r + ik e cot 9/r 2 + k 2 )]A' = (a + c )B' , 

(29) 

(iuj — ci + ic 2 k r + ic 3 kg + (3k 2 )B' 

= (ai — ia 2 k r — ia 3 kg + ak 2 )A' , (30) 

where k 2 = k 2 + kg/r 2 . The typical length scale in the radial direction is R & , while the 
typical angular scale in the latitude direction is tt, therefore we have k r = 2mn/R Q and kg = 
2im/ir = 2n, where m and n are non-zero integers. For the fundamental wave, m — n — 1. 
Consequently, k = 2[-K 2 m 2 + (R & n/r) 2 } 1 / 2 / R & . 

The dispersion equation can be obtained by multiplying Eq. (29) and (30). Defining 
Co = iuj + iU(l3/r + c 2 )k r + (f3/r 2 cot 9 + c 3 )k e ] 

—a + pk 2 , 

a 4 = a 2 k r + a 3 k e , 

c 4 = (l3/r + c 2 )k r + {[3 cot 9 / r 2 + c 3 ) k e , 

c 5 = (3(k r /r + kg cot 9/r 2 )(c 2 k r + c 3 kg), 

di = cl — cl + 4c 5 + 4(a + c )a 1: 

d 2 = Ac\(3(k r /r + kg cot 9/r 2 ) — 4(a + co)a4 — 2ciC4, 



we obtain 



Its solutions are 



Co 2 = ^(d 1 + id 2 ). (31) 



"2 



where 



D = (dj + ^)V2. (33) 
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The required solution is uj + because we need a positive real part of u to express the wave 
frequency: 

uj = uj + = u r + iui, (34) 

where 



W r = - 



(P/r + c 2 )k r - {(3 cot 6/r 2 + c 3 )k e 



00,; = 



(35) 
(36) 



The real part of a;, u r > 0, determines the period of the dynamo wave, and the imaginary 
part of ou, u>i, determines whether the dynamo wave grows (u>i < 0) or decays (u>i > 0). 



The period T is defined by 



and the growth rate T is defined by 



T = 27r/ccv, 



-OJi. 



These are two basic parameters of a dynamo. From their expressions we can see which 
mechanisms play a role. 



3.3. at-Vt dynamo 

When the buoyancy of the magnetic tube and the dynamic pressure of downflows are 
neglected, (5 = 0, the dispersion relation reduces to the usual version of the dynamo wave: 

u = ^y^EE- iy m^j, (37) 

where 

D = \a\ + (a 2 k r + ash) 2 } 1 ' 2 . (38) 

Since a±, a 2 and are related to Q, Eq. (37) shows that it is the a- and fi-effects that 
determine the period and growth rate of the usual dynamo. That is why one calls it a-VL 
dynamo. 

When Q is fixed by the observation, shown in Fig. 2, the other adjustable parameters are 
wavenumbers m and n, and the a parameter. Observations suggest that m — 1 and n = 2 
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for the equatorward branch, and m — — 1 and n = — 2 for the poleward branch, noticing that 
k r = 27im/R Q and k e = 2im/-n. The simplest a model is that a is a constant. The second 
simplest model is that a is some kind of function of radius and heliographic latitude. It is 
most likely that a depends upon dynamo magnetic fields. 



3.4. Anti-,3 dynamo 

The parameter (3 in the dynamo equations plays a role of dissipation. However, when 
j3 is negative, the dissipation effect reverses, we have an anti-dissipation effect, called the 
anti-/3 effect. When this effect dominates, the model is called an anti-/3 dynamo. 

In our model equations, we use (3 to explicitly include the buoyancy and downflows 
effects on a magnetic tube. When v z is negative and the reference field B is weaker than the 
critical field B c determined by fa + f3 v — 0, f3 is negative, we expect a positive growth rate 
for the tube field; when v z is negative and the reference field B is stronger than the critical 
field B c , (3 is positive, we expect that the tube field damps; when v z is positive, the tube 
field should always damp; if we ignore the downflows, the tube should move upwards and 
the field at the specified place dissipates quickly, depending how strong the tube field is. We 
want to examine if we have correctly modeled these intuitive results. 

The condition u r > suggests 

(P/r + c 2 )k r + (/3 cot 9/r 2 + c 3 )k e = 0. (39) 

If we demand ((3/r + c 2 )k r = and ((3 cot 9/r 2 + c 3 )ko = 0, we find that the reference field 
satisfies 

* " "5 ( 1 + !)"*'• (40) 

The amplitude is assumed to be proportional to the critical field B c . We can further calculate 
the second-order derivatives which appear in the expression of c\. 

^ = B[B 2 e + 1(1 +p/p b )csc 2 e + (P/p b )B e cot^], 
d 2 B B ( ~ 2 A 3 (3 V f. 1 



drl = 7>{ B '- B ' + ¥i B ' +1 ^- p ^ 

-i|(p-»(2 + P)-5|(2 + P) 
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1 L dp dP \ 1 A, ( dp dP 



8 \ <9 In r dlnr 8 3 h \ dlnr dlnr 



where p = din p/dlnr = r/H p and P = din P/ dlnr = r/H p . Here H p and H p are density 
and pressure scale heights respectively. 

We assume A = frB, where / is a constant factor. A negative / indicates that the 
poloidal field is in anti-phase with the toroidal field. Consequently, we obtain 

9A f Ft 



d6 

dA 
dr 

d 2 A „ d 2 B 



fB(l + B r ), 



fr- 



QQ2 ■> QQ2 ' 

d 2 A 2fB ^ | ^d 2 B 
dr 2 r dr 2 

These expressions are used to calculate c . 

Numerical calculations by using the standard solar model (Winnick et al. 2003) show 

that 

ci < 0, 

dx :=» c\ + Aa'ax > \d 2 \, 

d 2 w 4ci/?(A; r /r + k e cot 9/r 2 ) + 4a'a 4 , 

where a' = a + c . Since c depends upon B, a' depends upon B even if a does not depend 
upon B. The expression of Co given above provides us with a concrete nonlinear a model. 
Consequently, the general dispersion relation reduces to 

uj ~ | — P(k r /r + k e cot 6/r 2 ) + a'a A /ci\ 

+i((3k 2 + a'ai/ci + ^/16c?). (42) 



Letting a' = 0, T ps — /3/c 2 , which confirms the qualitative requirements stated above. 

In order to investigate the phase relationship between the poloidal and toroidal field, 
we set a = 0. Since a\ < near the surface, c\ < and the sign of Co is determined by the 
sign of /, we see that the term a'ai/ci > when / is positive and vice versa. This shows 
that the growth rate T ps —{/3k 2 + a'ai/ci) is larger when the two fields are in anti-phase 
with each other than when they are in phase. Obviously, this anti-phase relationship can be 
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attributed to the nonlinear a effect. Physically, it is due to the nonlinear interaction between 
the poloidal and toroidal fields. 

The observed equatorward and poleward branches of the solar magnetic field wave sug- 
gest 

f 2n(+l)/R Q if 9 > 45°, 
K = { 43 
I 2n(-l)/R G if < 45°, 

f 27r(+2)/7r if£>45°, 

kg = { 44 

{ 27r(-2)/7r if < 45°. 
The third term in the expression of the growth rate 

T = -((3k 2 + a'ai/ci + d 2 2 /l6cl) 
may explain the observed branching if we use the following a model: 



a 



+a if 9 < 45° and k r > and k e > 0, 

-ck if 9 < 45° and A; r < and k e < 0, 

-«i if 9 > 45° and k r > and > 0, 

+«i if ^ > 45° and k r < and < 0, 



(45) 



where ao and cti may slightly depend on 9, but may strongly depend upon r since the 
important ^-dependence of a has been written down explicitly. This is another nonlinear a 
effect, different from that specified by cq. In fact, a may be represented by the perturbation 
magnetic field B' as follows (e.g., Seehafer et al. 2003): 

1 



a = -r cor B'-(VxB'). 

Using the wave-like solution, Eq. (46), one can obtain Eq. (45). 

We have two options to qualitatively explain why the poleward branch is weaker that 
the equatorward branch: 

• «! < ao, 

• the observed poleward branch is the reflection wave, which is one part of the original 
wave. 



Because of the phase delay of the poleward branch with respect to the equatorward branch, 
the second option is preferable. 
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Solving Eq. (29) for A' Q , we obtain 



where 



Ri 
$1 



A' = Rie B', 

\a'\/[(d 2 2 /lQcl) 2 + a'a 4 / Cl ) 2 } 1/2 , 
arctan[a / a4/ci, — (d%/16cf)]. 



Using R\ and $i, we can express the temporal evolution of the toroidal magnetic field B 
and the poloidal vector potential A in terms of the initial field B for each magnetic flux 
tube as follows: 



B = B e 
A = B R ie Tt e m 



n exp[ ? (^-k-x+$ )] iffl>45°, 
exp [i(^ + k-x + $ )] if^<45°, 

exp [i -k-x + $ )] if^>45°, 
exp [i(^ + k-x + $ )j if#<45°, 



where <E> is determined by the initial condition. Substituting these expressions into Eq. (1), 
we obtain, 



B e 



rt 



where 



a 

a 4 



exp [i (M - k ■ x + $ )] [ e<j) + f e^a"] 

ae> 45°, 

exp [i + k • x + %)] [e + f e^a+] 

if < 45°, 



(cot 9 — ikg)e r — (1 — irk r )eg, 
(cot + ikg)e r — (1 + irk r )eg. 



(46) 



For the equatorward branch, we can investigate how its peak propagates with time by 
setting 

2nt/T - k r .r - k e 9 + % = (n - 1)tt, (47) 

where n equals to the integer part of 1 + 2t/T. Since k r = tt/Rq and kg = 4, we obtain the 
migration equation of the wave peak 



9 = 2tt/9 + [2t/T - (n - 1)]tt/4, 



(48) 



where we have assumed 9 = 40° when t — 0. Obviously, the integer n plays a role of cycle 
number. 
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For the poleward branch, we can only see the reflected wave because it propagates into 
the interior of the sun, k r = —2ti/R q . The reflection point is assumed to be located at the 
midpoint of the convection zone for the dynamo waves generated near the surface layers. This 
leads to r r fl m 0.86R Q . The phase velocity of the wave in the radial direction approximately 
equals to 100 cm/s. Therefore, it takes about to — 6 years for the wave to travel from the 
source region (r src = O.9995i? ) to the reflection point (r r fl = O.86-R ), and then to the 
surface (r = R Q ). As a result, the wave peak propagation equation for the poleward branch 
reads as follows: 

2n(t + t )/T + k r r + k e 9 + $ = (n' - 1)tt, (49) 
where n' equals to the integer part of l + 2(t + t )/T. The solution of this equation is 

= 5tt/18 - [2(t + t )/T - (n - 1)]tt/4, (50) 

where we have assumed 9 = 50° when t + t — 0. This implies a phase lag of ir/2 with 
respect to the equatorward branch, as shown in Fig. 3, which plots Eqs. (48) and (50). 



4. Comparison with observations 

In this section we try to compare the anti-/3 dynamo model with relevant observations. 
The model parameters are listed as follows: (i) the downflow velocity v z ~ —2.1 • 10 5 cm/s 
is fixed by the 3d numerical simulations (Robinson et al 2003); (ii) we adjust the flux tube 
radius a to obtain the observed cycle period T 22 years; (iii) £ = 1, which means that 
the flux tube is totally located within the downflow regions; (iv) the reflection point of the 
dynamo wave is assumed to be located at the midpoint of the convection zone; (v) the 
reference magnetic fields vary from zero to the critical field B c , which means that we should 
integrate all the observable quantities such as the cycle period over the reference field from 
zero to B c ; (vi) (3 t = \v 2 z (p/ ' Pfl 2 2a; (vii) a = 0; (viii) / = -lO" 10 . 



4.1. Period 



The /9-effect dominates in the anti-/3 dynamo model; the period T near the surface can 
be approximately expressed as follows 



T(a, 9) = 2tt/ 



' c -(3{2ti/R q + 4 cot 9/R 2 Q )dB/B c 

if 9 > 45°, 
J Bc {3(2tt/R g + 4 cot 9/R%)dB/B c 

if 9 < 45°. 
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If one fixes T « 22 years, one finds a=(1.15, 1.1, 1.05, ...,0.65, 0.59, 0.53, 0.45, 0.37, 0.27, 
0.15) Mm for 6>=(5°, 10°, 15°,..., 55°, 85°), respectively. This is comparable to 2 Mm, the 
typical size of flux tubes observed in the low latitude region. 



Similarly, the growth rate T of a flux tube near the surface can be approximately ex- 
pressed as follows 



where k 2 = (2-7T + 4) 2 /i?|,. For the parameters given above, one finds T 1 = (0.46, 0.48, 0.51, 
0.53, 0.56, 0.59, 0.62, 0.66, 0.71, 0.76, 0.82, 0.94, 1, 1.18, 1.44, 1.97, 3.54) years. This shows 
that the growth time of a typical flux tube of size a = 2 Mm is about several months. 

The flux tubes should stay in the convection zone long enough (i.e., for a time larger 
than or equal to their growth time) for them to be amplified by the dynamo action. Since 
our formulation has already included the buoyant dissipation effect, a positive growth rate 
measured by the growth time indicates that the downflows can inhibite the buoyant rise of 
the flux tube on the time scale specified by the growth time. Here we assume that the flux 
tube is so distorted that it is totally immersed in the downflow regions by setting £ = 1. 
From Fig. 5 of Fan et al (2003), we can see that only a small fraction of their flux tube 
is immersed in the downflow regions, which means that their £ parameter is much smaller 
than the value we use. Therefore, our results are not in conflict with their 3D simulation, 
which shows that the pumping has almost no effect in inhibiting the buoyant rise of the flux 
tube when the energy of the tube is larger than that of the surrounding environment. If the 
energy of the tube is comparable to that of the surrounding, it still rises to the surface, but 
being largely distorted. If we use a small £ parameter value, we can obtain similar results. 



Here we have two reasons to explain that the strong activity may be confined to low 
heliographic latitudes |$| < 35° of the equatorward branch, as observed: 

(i) growth rate at low latitude is larger than that at high latitude; 

(ii) what we observe at high heliographic latitudes are the reflected waves. They must be 
weaker than the original waves since the reflection is not necessarily complete. 



4.2. 



Growth rate 




4.3. Active regions 
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4.4. Phase dilemmas 

(1) Eq. (46) shows that the phase difference between the poloidal and toroidal field is $i. 

Calculations show $ x m 170° for 0° < 9 < 180°. This is in good agreement with the 
observation: the poloidal field is almost in antiphase with the toroidal component; 

(2) The 7r/2 phase lag of the poleward branch with respect to the equatorward branch can 

be explained by the time lag of 6 years for the former to travel from the source region 
to the reflection point, and then from the reflection point to the surface. 

4.5. Butterfly diagram 

Fig. 3 depicts the latitude of the peak of the dynamo wave as a function of time. Only 
when we select the mid-point of the convection zone as the reflection point, the time delay 
of the poleward branch is about 6 years. This is a bit arbitrary. 

4.6. Cyclic variations of solar global parameters and helio-seismic frequencies 

We construct models of the structure and evolution of the Sun which include variable 
magnetic fields and turbulence. The magnetic effects are: (1) magnetic pressure, (2) mag- 
netic energy, and (3) magnetic modulation to turbulence. The effects of turbulence are: (1) 
turbulent pressure, (2) turbulent kinetic energy, and (3) turbulent inhibition of the radiative 
energy loss of a convective eddy, and (4) turbulent generation of magnetic fields. Using these 
ingredients we construct five types of solar variability models (including the standard solar 
model) with magnetic effects. These models are in part based on three-dimensional numerical 
simulations of the superadiabatic layers near the surface of the Sun. The models are tested 
with several sets of observational data, namely, the changes of (1) the total solar irradiance, 
(2) the photospheric temperature, (3) radius, (4) the position of the convection zone base, (5) 
low- and medium-degree solar oscillation frequencies. We find that turbulence plays a major 
role in solar variability and only a model which includes a magnetically-modulated turbulent 
mechanism can agree with all the current available observational data. We find that because 
of the somewhat poor quality of all observations (other than the helioseismological ones), 
we need all data sets in order to restrict the range of models. These studies, presented in 
Li et al (2003), strongly suggest a near-surface dynamo model. The anti-/3 dynamo, which 
operates near the surface, may meet this need. 
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5. Conclusions 

Helioseimic observations of the cyclic variations of low- and medium-degree solar os- 
cillation frequencies are very precise. The cyclic variation of the position of the convection 
zone base, obtained by helioseismic inversions (Basu & Antia 2000; Basu 2002), is also very 
accurate. These high-precision observations seem to exclude the location of strong magnetic 
fields at the base of the convection zone or in the overshoot layer, as is generally assumed on 
the basis of stability against magnetic buoyancy and equipartition arguments. This is one 
of the main motivations to reconsider dynamos located in the convection zone. The second 
motivation comes from trying to reproduce the observed cyclic variations of the total solar ir- 
radiance (Frohlich & Lean 1998), the solar effective temperature (Gray & Livingston 1997), 
and the solar radius (See §2.3 of Li et al 2003 for a brief review). Extensive numerical ex- 
periments (Li & Sofia 2001) show that only magnetic fields located in the convection zone 
can simultaneously produce the observed cyclic variations of these three global solar param- 
eters. The third motivation originates from numerical simulations of the solar convection 
zone (Chan and Sofia 1989; Nordlund 1999; Robinson et al 2003), because they indicate a 
major presence of downward-moving plumes of high velocity. 

However, convection zone dynamos have to face the magnetic buoyancy instability prob- 
lem. The best way to solve this problem is to incorporate the buoyancy into the dynamo 
equations. This forces us to treat the mean field in a flux tube. To include the dynamic 
pressure of a flow, the flux tube should be regarded as a body. As a field, both the cyclonic 
convection (a-effect) and differential rotation (fi-effect) play a role. As a body, the tube ex- 
periences not only a buoyant force, but also a dynamic pressure of downflows above the tube. 
We show that these two dynamic effects can be incorporated into the dynamo equations by 
converting them in two diffusion terms. 

We analyzed and solved the extended dynamo equations in the linear approximation by 
adopting the observed solar internal rotation rate, and by assuming a downflow suggested 
by numerical simulations of the solar convection zone. The results are in agreement with the 
basic properties of the solar cycle, and give a theoretical basis to understand the variations 
of solar irradiance, temperature, radius and helioseismic frequencies in the course of the 
magnetic cycle. 

Finally, the model results give a theoretical basis to interpret the observed cyclic varia- 
tions of the total solar irradiance, the solar effective temperature, the solar radius, and the 
helioseismic frequencies (Li et al. 2003; Le et al 2002; Li & Sofia, 2001; Sofia & Li, 2001). 

We want to thank Dr. Sarbani Basu for her kindly providing us with rotation rate data, 
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A. A. Dynamo equations in spherical coordinates 



In this paper we assume azimuthal symmtry, which means that all quantities do not 
depend upon the azimuthal coordinate 0. 



Using Eq. (1) to expand Eq. (16), we obtain 



dB 
~dt 



■e^+Vx 



- «B - (3V 2 A ) e 4 



Vx{^x [Vx (Ae^+aVx (Ae^-pVx (Be*)}, 

(Al) 



where we have used 



U = ft x r = rfl sin 9e^ = Ue<h. 



(A2) 



This preprint was prepared with the AAS IATeX macros v5.2. 
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We have also used the mathematical equality 

VxVx(^) = V(V • Ae^) - V 2 (Ae^) 

= -e V 2 A (A3) 

The second equality hold well when A does not depend upon <p. 

We use the following mathematical equality to calculate the curl of a vector field F, 



VxF 



r sin 9 



9 c ^ \ dF <> 



e r +- 



sin# d<f) dr {r ^ 



e e +- 



d dF r 
d~r {rFe) ~ W 



Since we assume that the system is in azimuthal symmetry, we know 

V x (Fe ) = -_L.A( s in0F)e r - -|-(rF)e e 
r sin 9 39 r or 

for F = F(r,9)e ( j ) . Using these formula in Eq. (Al), we obtain 



dB M\ 



dA 



f3V 2 A e, 



0, 



(A4) 
(A5) 

(A6) 



where 
M = 



rV x 
+ 



U 9 , a d . . „ (3 d 

-jr( rA ) + — smM 

r or r sin 9 o9 



U 



d ad 3d 

(sinM) - -j-irA) + '-j-irB) 
r sin 9 09 r Or r Or 



r sin 9 89 



sin 9B) 



d_ 

dr 



JL 9 {sin9A ) - a^{rA) + (3^{rB) 
sin 9 o9 or or 



d_ 
89 



r ar r sin tl« r sin y at/ 



sin 05) 



dr 



dA dA dB 

U— + (U cot - a) A - ar— + f3B + f3r— 
d9 dr dr 



d_ 

m 

d_ 

dr 



ydA adA 
dr r d9 



U + acot9 . (3dB (3cot9„ 

A - tttt - B 

r r d9 r 



(17 cot 0- a) - -^-[U + acote)] }A+ I (U cot 6- 
r d9 J [ 



a 



d~r iar) - M 



d\J\ dA 

dr 



U + a cot 9 d fa 
r d9 1 Vr 



dU] dA 



dr 



d 2 A a d 2 A 



+ 



d{3 
dr 



d_ 





B + 







d9 

d 



ar 



dr 2 
dB 



dr 



+ 



r d9 2 
(3 cot 6 



dB 



n d 2 B (3d 2 B 
jf.p r 1- _ . 

dr 2 r d9 2 



(AT) 
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Eq. (A5) shows that the second term of Eq. (A6) on the left hand side is perpendicular to 
the first term. Consequently, Eq. (A6) is equivalent to the following two scalar equations: 



OA 
~dt 



aB + p 



1 dA cot 9 dA d 2 A 1 d 2 A 



r dr 



r 2 d9 dr 2 



de 2 



dB . OA dA 



a 



d 2 A 
dr 2 



+ 



d 2 B 1 8 2 B 
+ 



dr 2 



de 2 



(A8) 



1 d 2 A\ 



, „ , dB , dB 
+ b 1 B + b 2 — + b 3 — 



dr 



89 



(A9) 



where we have used 



V 2 A 



l_d_ 

r 2 dr 
Id A 
r dr 



,8A 
dr 



+ 



d 



r 2 sin 9 d9 



sin# 



dA 

~d~e 



cotd d A d 2 A Id 2 A 



r 2 89 dr 2 



de 2 



(A10) 



and defined 



ai 



a 2 = 



a 3 = 



b 2 
bs 



d 



1 d 



, (Ucot9-a)--—(U + acot9) 
dr r 89 

n a ^a ■ n a 9/i da 1 ^dcA 

iZcos 6AZ r — ililg sm9 H — esc — cot 9— 

r dr r 86 J 



(UcotO- a) 
U + a cot 6 



8_ 

dr 
d 



ar 



dU_ 



d A + _ 
dr dd 



r + d9 
d ( (3 cot 9 



V r / <9r J r 



= — fif^ sin 

fifL sin 9 — 



2a da 
r dr ' 
a cot 6* 



2pb (dB_ cot 9 dB 
rB \ dr r d9 



r 

,2 , 



1 da 
rd9 



,(A11) 



pcsc 2 9 p 



r B dr r ' 



P cot 



— [pcot9 + 



2p b dB 
~B~89 



The valid condition for the second equality sign is B <C \^8nP. We define fi r = dlnQ/dlnr 
and fig = d hi Q/ 86. We also define 



P 



' dp 
d\nr 



Pt (dhip _ d\nP 
2 \ d In r <9 In r 



-A 



3<91nP l<9mp 

2 H - 

2 <91nr 



A, 

2<91nr/ 2 



<91nP <91np s 
<91nr <91nr 
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B. B. Linearization of Dynamo equations 

In order to linearize Eqs. (17-18), we decompose B and A into two parts: 

B = B + B', (Bl) 
A = A + A', (B2) 

where B and A represent the reference state of the system, while B' and A' represents the 
perturbations of B and A near Bq and Aq. Substituting Eqs. (B1-B2) into Eqs. (17-18), we 
obtain two linear equations that govern B' and A', 

OA' . n fldA' cot 9 d A' d 2 A' 1 d 2 A'\ 

(a + c )B> + &-_ + __ + _ + __ , (B3) 



dt \r dr r 2 89 dr 2 r 2 d9 2 

OB' dA' a dA' (d 2 A 1 d 2 A'\ 

~~dT = a ' A+a2 ^ + a "W~ a {-d^ + 7 2 ^9 2 -) 
dB' dB' a fd 2 B' i d 2 B'\ 

+ ° 2 ~dV + Cs ~d9~ + ydr^ + ^ 2 ~d9 2 ~ ) ' (B4) 

where the superscript or the subscript o means that the quantities are evaluated at the 
reference state, for example, a\ = ai\ b=b ,a=a , and 



c 



2$ fldA cot9dA d 2 A 1 d 2 A 
B \r dr r 2 89 dr 2 r 2 d9 2 
,o dh °i o db 2 dB o db° 3 dB 2/3° fd 2 B 1 d 2 B { 



a = Vt + —±-B 



dB u dB dr 8B 89 B \ dr 2 r 2 89 

db\ db° 2 dB , <%° dB 

c 2 = o 2 + ——B 



dB' u dB' dr dB' 89 ' 
7() db\ , dbl dB Q , «%° dB 



C 3 = b 3 + TT^ B 0+ ao „ o + 



<9B ' u <95 ' dr dB'> 89 



Here _B = dB /dr and B ' = dB /d9. Using the expressions for b's given in Appendix A, 
we can calculate their derivatives with respect to Bo, B' , and Bq\ 

dbl _ 2/3° (dBp | cot OdB \ 2$ esc 2 9 2{3 b ( | 3<91nP ldlnp 



dB rB 2 \ dr r d9 J B r 2 r 2 B \ 2d In r 2<91nr 
<%°, _ 4/3° , 2 / 3°<9B 2/3 6 / 2+ 3 91nP ldlnp 



<9£> ^-Bo -^o ^ r rB o V 2(9 In r 2 9 In r 



dB r 2 \ B B 2 d9 
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dB'o " rBo 1 



dB' 
d% 
dB' 
db\ 
dB» 

dB>> 

W_ 
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2^ 
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2$ cot e 
tBq r 

0, 

r 2 Bn 




Fig. 1. — Schematic of the force balance of a magnetic tube in the convection zone. 
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Fig. 3. — Heliographic latitude of the dynamo wave peak as a function of time 



